STEP 5: Plot

Plotting multispectral data

Multispectral data can be plotted as:

  1. Individual bands
  2. Spectral indices
  3. True color 3-band images
  4. False color 3-band images

Spectral indices and false color images can both be used to enhance images to clearly show things that might be hidden from a true color image, such as vegetation health.

Try It: Import libraries

Add missing libraries to the imports

import cartopy.crs as ccrs # CRSs
# Interactive tabular and vector data
import hvplot.xarray # Interactive raster
# Overlay plots
import numpy as np # Adjust images
import xarray as xr # Adjust images
See our solution! 
import cartopy.crs as ccrs # CRSs
import hvplot.pandas # Interactive tabular and vector data
import hvplot.xarray # Interactive raster
import matplotlib.pyplot as plt # Overlay plots
import numpy as np # Adjust images
import xarray as xr # Adjust images

There are many different ways to represent geospatial coordinates, either spherically or on a flat map. These different systems are called Coordinate Reference Systems.

Try It: Prepare to plot

To make interactive geospatial plots, at the moment we need everything to be in the Mercator CRS.

  1. Reproject your area of interest with .to_crs(ccrs.Mercator())
  2. Reproject your NDVI and band raster data using `.rio.reproject(ccrs.Mercator())
See our solution! 
# Make sure the CRSs match
aoi_plot_gdf = denver_redlining_gdf.to_crs(ccrs.Mercator())
ndvi_plot_da = denver_ndvi_da.rio.reproject(ccrs.Mercator())
band_plot_dict = {
    band_name: da.rio.reproject(ccrs.Mercator())
    for band_name, da in band_dict.items()
}
ndvi_plot_da.plot(cmap='Greens', robust=True)
ndvi_plot_da.hvplot(geo=True, cmap='Greens', robust=True)

Plot raster with overlay

Try It: Plot raster with overlay using xarray and geopandas

Plotting raster and vector data together using pandas and xarray requires the matplotlib.pyplot library to access some plot layour tools. Using the code below as a starting point, you can play around with adding:

  1. Labels and titles
  2. Different colors with cmap and edgecolor
  3. Different line thickness with line_width

See if you can also figure out what vmin, robust, and the .set() methods do.

ndvi_plot_da.plot(vmin=0, robust=True)
aoi_plot_gdf.plot(ax=plt.gca(), color='none')
plt.gca().set(
    xlabel='', ylabel='', xticks=[], yticks=[]
)
plt.show()
See our solution! 
ndvi_plot_da.plot(
    cmap='Greens', vmin=0, robust=True)
aoi_plot_gdf.plot(
    ax=plt.gca(), 
    edgecolor='gold', color='none', linewidth=1.5)
plt.gca().set(
    title='Denver NDVI July 2023',
    xlabel='', ylabel='', xticks=[], yticks=[]
)
plt.show()

Try It: Plot raster with overlay with hvplot

Now, do the same with hvplot. Note that some parameter names are the same and some are different. Do you notice any physical lines in the NDVI data that line up with the redlining boundaries?

(
    ndvi_plot_da.hvplot(
        geo=True,
        xaxis=None, yaxis=None
    )
    * aoi_plot_gdf.hvplot(
        geo=True, crs=ccrs.Mercator(),
        fill_color=None)
)
See our solution! 
(
    ndvi_plot_da.hvplot(
        geo=True, robust=True, cmap='Greens', 
        title='Denver NDVI July 2023',
        xaxis=None, yaxis=None
    )
    * aoi_plot_gdf.hvplot(
        geo=True, crs=ccrs.Mercator(),
        line_color='darkorange', line_width=2, fill_color=None)
)

Plotting bands as subplots

Try It: Plot bands with linked subplots

The following code will make a three panel plot with Red, NIR, and Green bands. Why do you think we aren’t using the green band to look at vegetation?

raster_kwargs = dict(
    geo=True, robust=True, 
    xaxis=None, yaxis=None
)
(
    (
        band_plot_dict['red'].hvplot(
            cmap='Reds', title='Red Reflectance', **raster_kwargs)
        + band_plot_dict['nir'].hvplot(
            cmap='Greys', title='NIR Reflectance', **raster_kwargs)
        + band_plot_dict['green'].hvplot(
            cmap='Greens', title='Green Reflectance', **raster_kwargs)
    )
    * aoi_plot_gdf.hvplot(
        geo=True, crs=ccrs.Mercator(),
        fill_color=None)
)
See our solution! 
raster_kwargs = dict(
    geo=True, robust=True, 
    xaxis=None, yaxis=None
)
(
    (
        band_plot_dict['red'].hvplot(
            cmap='Reds', title='Red Reflectance', **raster_kwargs)
        + band_plot_dict['nir'].hvplot(
            cmap='Greys', title='NIR Reflectance', **raster_kwargs)
        + band_plot_dict['green'].hvplot(
            cmap='Greens', title='Green Reflectance', **raster_kwargs)
    )
    * aoi_plot_gdf.hvplot(
        geo=True, crs=ccrs.Mercator(),
        fill_color=None)
)

Color images

Try It: Plot RBG

The following code will plot an RGB image using both matplotlib and hvplot. It also performs an action called “Contrast stretching”, and brightens the image.

  1. Read through the stretch_rgb function, and fill out the docstring with the rest of the parameters and your own descriptions. You can ask ChatGPT or another LLM to help you read the code if needed! Please use the numpy style of docstrings
  2. Adjust the low, high, and brighten numbers until you are satisfied with the image. You can also ask ChatGPT to help you figure out what adjustments to make by describing or uploading an image.
rgb_da = (
    xr.concat(
        [
            band_plot_dict['red'],
            band_plot_dict['green'],
            band_plot_dict['blue']
        ],
        dim='rgb')
)

def stretch_rgb(rgb_da, low, high, brighten):
    """
    Short description

    Long description...

    Parameters
    ----------
    rgb_da: array-like
      ...
    param2: ...
      ...

    Returns
    -------
    rgb_da: array-like
      ...
    """
    p_low, p_high = np.nanpercentile(rgb_da, (low, high))
    rgb_da = (rgb_da - p_low)  / (p_high - p_low) + brighten
    rgb_da = rgb_da.clip(0, 1)
    return rgb_da

rgb_da = stretch_rgb(rgb_da, 1, 99, .01)

rgb_da.plot.imshow(rgb='rgb')
rgb_da.hvplot.rgb(geo=True, x='x', y='y', bands='rgb')
See our solution! 
rgb_da = (
    xr.concat(
        [
            band_plot_dict['red'],
            band_plot_dict['green'],
            band_plot_dict['blue']
        ],
        dim='rgb')
)

def stretch_rgb(rgb_da, low, high, brighten):
    """ 
    Contrast stretching on an image

    Parameters
    ----------
    rgb_da: array-like
      The three channels concatenated into a single array
    low: int
      The low-end percentile to crop at
    high: int
      The high-end percentile to crop at
    brighen: float
      Additional value to brighten the image by

    Returns:
    --------
    rgb_da: array-like
      The stretched and clipped image
    """
    p_low, p_high = np.nanpercentile(rgb_da, (low, high))
    rgb_da = (rgb_da - p_low)  / (p_high - p_low) + brighten
    rgb_da = rgb_da.clip(0, 1)
    return rgb_da

rgb_da = stretch_rgb(rgb_da, 2, 95, .15)
rgb_da.plot.imshow(rgb='rgb')
rgb_da.hvplot.rgb(geo=True, x='x', y='y', bands='rgb')

False color images

Try It: Plot CIR

Now, plot a false color RGB image. This is an RGB image, but with different bands represented as R, G, and B in the image. Color-InfraRed (CIR) images are used to look at vegetation health, and have the following bands:

  • red becomes NIR
  • green becomes red
  • blue becomes green
Looking for an Extra Challenge?: Adjust the levels

You may notice that the NIR band in this image is very bright. Can you adjust it so it is balanced more effectively by the other bands?

See our solution! 
rgb_da = (
    xr.concat(
        [
            band_plot_dict['nir'] * .8,
            band_plot_dict['red'],
            band_plot_dict['green']
        ],
        dim='rgb')
)

rgb_da = stretch_rgb(rgb_da, 2, 98, 0)
rgb_da.plot.imshow(rgb='rgb')
rgb_da.hvplot.rgb(geo=True, x='x', y='y', bands='rgb')